Testing passed. I welcome any thoughts or review from @paciorek or @perrydv

So do we think the most common use case would be that y[1:n] is a leaf node containing some data and some non-data elements and therefore the non-data are predictive nodes? In that case it looks like we assign a RW_block sampler. So @danielturek I'm not sure what you meant by "Any latent dimensions which have no stochastic dependencies are assigned posterior_predictive samplers." More substantively, should we work out the true MVN conditional distribution in that situation and sample from the true predictive? If we do think this is a common case, then I think it might be worth it. (And I would be happy to lend a hand if we decide to do it.)

@paciorek Please see the addition of the partial_mvn_pp sampler ("partial MVN posterior predictive" sampler), which is used internally in the partial_mvn sampler.

When one or more dimensions (from among the unobserved dimensions) of a partially-observed MVN node themselves have no data dependencies, then a partial_mvn_pp sampler is now assigned to these dimension(s) - instead of either RW or RW_block samplers, as before. This new partial_mvn_pp generates draws directly from the conditional MVN distribution.

@paciorek A second set of eyes would be appreciated. There might also be some efficiency improvements with the linear algebra within the partial_mvn_pp sampler, which could be improved upon from the present implementation.

@paciorek Thanks for the careful review.

Repeating the request from above, my first preference would be to make the backsolve and forwardsolve changes myself. But that would take a bit more time for me to understand and make sure I do it right, and I'm out of time for today. If you are willing to make this change yourself, and leave the original (current) code there commented out, I would very much appreciate that, and I'll take a very careful look later.

I'd be happy to discuss strategies for "caching" the components of the covariance matrix. But I think it would have to be done with some thought and some caution. Even if the covariance node has no parent nodes (stochastic or otherwise), the contents of the node can still be modified (regardless if the value is provided by constants, data, or initial values). So I can imagine a few ways we might try to do this, for example:

"If the covariance node has no parent nodes, then we have a indicator internal to the sampler (say firstRun) which indicates if this is the first run of the MCMC (initially TRUE), and when firstRun is TRUE, we pull the value of the covariance from the model, and perform the factorizations etc, and cache those values in the sampler. Then, on subsequent runs of the sampler, when firstRun is FALSE, we just used the cached values. However, when the reset method is called, indicating starting a new run of the MCMC, we set firstRun back to TRUE, in case the user changed the value of the covariance (inside the model) in between MCMC runs. (and, we also take it on blind faith that a user didn't change the value of the covariance, then continue an MCMC run using reset = FALSE, which would make the entire MCMC a little questionable, and make the subsequent MCMC results incorrect in some sense)."

This is a tricky one, unless you have some better ideas that I'm not imagining.

In the case of assigning the barker sampler in place of the RW_block sampler, need we do anything to enable derivatives in the model?

@paciorek Yes, please go ahead and make / add any tests for any of the following:

• any testing with changing the nimble option MCMCuseBarkerAsDefaultMV to TRUE, and correct assignment and operation of the barker sampler.
• check if partial_mvn_pp is correctly assigned, build the MCMC in some of the tests and check the internal sampler assignments
• test that checks correctness of the partial_mvn_pp sampler in a basic case.

Ok I pushed a commit that has the implementation of the Cholesky of Sigma22. I didn't run tests (though CI will of course) or check correctness yet.

As far as Barker, yes, that's correct the model would need to have derivs built. But I don't think anything is different with using Barker here than in our default MCMC config if that Barker option is set.

Regarding caching, I will comment on that in a moment below.

I will try to work on some of the other things this afternoon as well. Or tomorrow.

Unless Frederic indicates urgency in response to my email I am now planning that we include this in 1.4.1 and am thinking I will push this to CRAN on Tuesday. I will be out of town over the long weekend.

Regarding caching, I think the case that we can safely use cached values is if the covariance is RHS-only. If it has no parent nodes it could still be being sampled in an MCMC (e.g., if it is dwishart or dlkj with fixed prior parameters).

Two other thoughts:

• we could do similar machinations with the mean so as to save recalculation of the conditional mean, which involves matrix-vector manipulations. Clearly less important but if we are working out the caching perhaps we might as well.
• we could in cases where we don't know if the values might have changed, actually check if the cached covariance (and perhaps mean) are the same as the current values from getParam. Then if they are, we can avoid some calculations. Checking should be O(n^2), while the saved calculations are O(n^3) (dominated by the Cholesky of either Sigma11 or Sigma22). But in cases where the covariance is changing every iteration then we are doing a bunch of useless checking. So I think it probably doesn't make sense to do this.

I guess if my reasoning about checking for RHS-only is correct, I propose we combine that with your sketch of using firstRun. And I think we might as well check if the mean is RHS-only too.

I'm working on some additional tests.

@danielturek I don't understand the logic in "Values change at same rows for unobserved parts of partially observed node when multivariateNodesAsScalars = FALSE" test. Based on the model is uses the _pp sampler, so multivariateNodesAsScalars is not involved. Don't you want a model where y is intermediate?

@paciorek That test predated the addition of the partial_mvn_pp sampler - meaning, that test existed before today. The idea behind that test was to confirm that all unobserved dimensions of the MVN node change in unison (or, none of them change) on each MCMC iteration, since when multivariateNodesAsScalars = FALSE, the unobserved dimensions should have a RW_block sampler assigned (and hence all change simultaneously, or all remain unchanged, on each MCMC iteration). This is of course now moot, since for the model in this test the partial_mvn_pp now updates all unobserved dimensions on every MCMC iteration, so they will all (always) change on each iteration. I'm not wedded to this test, there's more of a backstory here than I'll get into, and I would be happy to have this test be removed, or modified, if you're so inclined.

@paciorek I worked through the linear algebra. I imagine it's simple to you, but this was helpful for me, and it was nice to see especially how the forwardsolve can be used to calculate the conditional posterior for Sigma, leveraging the transpose-relationship between Sigma12 and Sigma21. I now have a much better understanding of the uses of forwardsolve and backsolve. I also fixed an extremely minor dimension mistake.

Ok, so my new code had some issues with dimensions that I just fixed. Results are now the same as Daniel's implementation.

I added some tests, but the one I still want to add is one to check correctness of the _pp samples. I'll do that tomorrow morning.

I'll let you decide what to do about the 'moot' test.

I think we stepped on each other's commits a bit, but I think everything is ok now.

@paciorek expect_success does not operate as one might imagine. Use expect_no_error instead.